Abstract

Markovian pure jump processes model a wide range of
phenomena, including chemical reactions at the molecular level, dynamics
of wireless communication networks and the spread of epidemic diseases
in small populations. There exist algorithms like Gillespie's SSA or
Anderson's Modified Next Reaction Method, that simulates a single
trajectory with the exact distribution of the process, but this can be
time consuming when many reactions take place during a short time
interval. Gillespie's approximated tau-leap method, on the other hand,
can be used to reduce computational time, but it may lead to
non-physical values due to a positive one-step exit probability, and it
also introduces a time discretization error. Here, we present a novel
hybrid algorithm for simulating individual trajectories which adaptively
switches between the SSA and the tau-leap method. The switching
strategy is based on a comparison of the expected inter-arrival time of
the SSA and an adaptive time step derived from a Chernoff-type bound for
the one-step exit probability. Because this bound is non-asymptotic,
we do not need to make any distributional approximation for the tau-leap
increments. This hybrid method allows us (i) to control the global exit
probability of any simulated trajectory, and (ii) to obtain accurate
and computable estimates of the expected value of any smooth observable
of the process with minimal computational work. We present numerical
examples that illustrate the performance of the proposed method.